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Abstract-  A  method  for  baseline  (BL)  removal  in  needle  EMG 
records  is  presented.  Different  processing  techniques  are 
sequentially  used.  Firstly  motor  unit  action  potentials  (MUAPs) 
are  extracted  from  the  signal  by  means  of  a  wavelet  transform- 
based  procedure.  Potential-free,  discontinuous  segments  are  thus 
obtained,  whose  fluctuation  is  assumed  to  be  related  to  BL 
wander.  These  signals  are  then  time  averaged  to  attenuate  the 
effect  of  noise  and  of  low  amplitude  MUAPs  originated  distant 
from  the  electrode.  Spline  interpolation  is  then  used  to  build  a 
continuous  reconstructed  signal  whose  spectral  characteristics 
approximate  that  of  the  real  BL.  The  spectrum  of  this  signal  is 
estimated  by  AR  modeling  and  an  FIR  filter  is  implemented 
accordingly  for  filtering  out  the  BL  low  frequency  components 
from  the  original  EMG  signal.  Two  merit  figures  are  devised, 
which  measure  the  degree  of  BL  fluctuation  present  in  an  EMG 
record.  These  figures  are  used  to  compare  our  method  with  the 
conventional  approach  which  consider  the  BL  to  be  a  constant 
value.  Experiments  for  BL  removal  from  real  and  simulated 
EMG  signals  are  carried  out.  The  superior  performance  of  our 
approach  is  shown  regarding  these  merit  figures  and  visual 
inspection. 
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I.  Introduction 

Motor  unit  action  potential  (MUAP)  expresses  the 
electrical  activity  of  the  muscle  fibers  of  a  motor  unit  (MU) 
recorded  from  a  needle  electrode.  The  shape  of  MUAP 
waveforms,  and  their  similarity  in  consecutive  appearances 
contain  valuable  information  about  the  state  of  a  muscle, 
helping  to  distinguish  between  myopathic,  neurogenic  or 
normal  states  and  to  measure  the  degree  of  abnormality. 
MUAP  analysis  is  thus  a  daily-work  procedure  in  clinical 
electromyography  (EMG).  For  MUAP  characterization 
different  parameters  are  used:  duration,  amplitude,  area, 
number  of  turns  and  phases,  jiggle,  etc.  [1],  [2].  The 
measurement  of  these  parameters  is  influenced  by  baseline 
(BL)  fluctuation  along  the  recorded  signal.  Particularly 
relevant  is  this  influence  on  duration,  and  jiggle  (as  measured 
by  the  CAD  parameter  [2],  [3]),  as  these  include  in  their 
definition  amplitude  criteria  with  respect  to  the  BL.  Therefore 
a  precise  estimation  of  the  BL  will  help  making 
measurements  of  these  parameters  more  accurate. 

Current  methods  consider  the  BL  to  be  constant  throughout 
the  MUAP:  they  either  average  the  samples  in  the  segments  at 
both  ends  of  the  MUAP  discharge  (typically  3  ms  segments 
are  taken  in  25  ms  registers  in  which  the  main  spike  of  the 
MUAP  occupies  the  central  position)  [1];  or  else  just  give  the 
BL  a  zero  level  (system  ground). 

We  regard  the  BL  as  a  low  frequency  fluctuation  present  in 
the  recorded  signal  (Fig.  1)  due  to  artifacts  of  different  nature: 
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the  movement  of  the  recording  needle  relative  to  the  muscle, 
variation  of  skin  potential  induced  by  the  needle,  electrical 
drifts  in  the  acquisition  equipment,  etc.  But  above  all,  the 
main  source  of  BL  fluctuation  is  the  activity  of  distant  MUs, 
which  produce  potentials  that  cannot  even  be  identified  as 
MUAPs,  as  they  only  provoke  a  mild  BL  wander.  In  the  EMG 
signal,  there  also  appear  “secondary”  potentials,  not  to  be 
confused  with  the  previous.  They  originate  relatively  far  from 
the  electrode  and  are  low-amplitude  and  smooth,  and  thus  not 
valid  for  EMG  clinical  analysis,  but  still  considered  as  real 
MUAPs  and  not  as  BL  components. 

Several  methods  for  BL  removal  have  been  applied  to  other 
biomedical  signals,  such  as  ECG  [4],  [5].  They  heavily  rely 
upon  “a-priori”  knowledge  of  the  BL  frequency  band.  In 
EMG,  BL  frequency  components  are  not  too  well  known  and 
are  subject  to  high  variability  across  different  muscles, 
individuals  and  recordings.  In  the  method  we  present  in  this 
paper  several  signal  processing  techniques  are  used  to 
estimate  the  spectrum  of  the  BL  present  in  an  EMG 
recording.  A  convenient  filter  is  then  used  to  filter  out  the  BL 
frequency  components. 

II.  Material 

Twenty  EMG  real  signals  from  the  right  tibialis  anterior 
muscle  of  a  39  year-old  healthy  man  were  analyzed.  The 
EMG  signals  were  recorded  at  different  degrees  of  voluntary 
contraction  using  an  electromyograph  (Counterpoint,  Dantec 
Co.,  Denmark)  and  disposable  concentric  needle  electrodes. 
After  antialiasing  filtering,  signals  were  sampled  at  2.4  kHz  (a 
higher  sampling  frequency  was  not  required  in  this  study). 

Ten  simulated  signals  were  also  analyzed.  Real  MUAP 
waveforms  and  secondary  potentials  were  taken  as  templates 
to  form  MUAP  trains.  Templates  were  repeated  at  a 
determined  frequency  (between  3  and  12  Hz)  for  each 


Figure  1 .  EMG  recording:  several  MUAPs  appear  on  top  of  a 
time-varying  baseline  contamination. 
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different  train)  which  was  not  constant  but  subject  to  a  small 
random  variation  (±2  Hz).  Simulated  signals  were  composed 
of  several  MUAP  trains  (between  1  and  5)  and  corrupted  by  a 
secondary  potential,  white  gaussian  noise  (SNR  between  30 
and  43  dB)  and  BL  fluctuations  (SNR  between  14.6  and  28.7 
dB).  The  BL  was  simulated  by  filtering  zero-mean  white 
gaussian  signals  with  a  low-pass  Butterwoth  filter  [6]  (cut-off 
frequencies  between  5  and  20  Hz).  The  characteristic  of  both 
signal  sets  exhibited  wide  variation  with  respect  to  the  degree 
of  activity  (number  of  present  MUAPs),  MUAP  duration  and 
spectral  range  of  BL  fluctuation. 

III.  METHODS 

The  method  comprises  several  sequential  phases: 

a)  estimation  of  EMG  activity  level; 

b)  segmentation  of  the  EMG  signal  into  segments 
containing  MUAPs  and  free-MUAP  segments  (BL); 

c)  BL  spectral  characterization; 

d)  BL  filtering. 

In  each  phase,  different  alternatives  were  tested  and  those 
which  yielded  better  results  were  selected.  In  the  sequel  we 
describe  these  phases. 

A.  Estimation  of  the  level  of  EMG  activity 

The  number  of  MUAPs  present  in  an  EMG  recording  (i.e., 
level  of  EMG  activity)  depends  on  the  contraction  level  of  the 
muscle.  With  a  low  contraction  level,  few  MUs  are  activated 
and  the  record  presents  relatively  long  and  steady  BL 
segments,  with  scattered  MUAPs.  As  the  contraction  force 
increases,  more  MUs  are  recruited  and  the  level  of  EMG 
activity  increases  while  the  BL  segments  without  MUAPs 
become  fewer  and  shorter.  We  use  the  number  of  spikes  per 
ms  as  a  quantitative  measure  of  the  level  of  EMG  activity. 

B.  MUAP  detection 

This  phase  focuses  on  the  MUAP  detection  and  isolation 
from  BL  segments,  making  use  of  the  discrete  wavelet 
transform  (DWT).  Specifically  the  non-orthogonal  quadratic 
spline  wavelet  is  used,  which  has  been  successfully  utilized 
for  the  detection  of  characteristic  points  in  ECG  signals  [7]. 
The  detection  algorithm  is  similar  to  that  presented  in  [7]  and 
its  applicability  is  also  sustained  by  the  fact  that  the  uniphasic 
(only-one -peak)  shape  of  this  wavelet  resembles  that  of  the 
uniphasic  basic  components  conforming  the  EMG  waveform. 
Similarly  to  the  QRS  detection  in  EEG,  MUAP  maxima, 
minima  and  zero-crossing  points  are  detected  in  the  EMG 
signal,  and  from  these,  MUAP  initial  and  end  points.  This 
part  of  the  process  is  further  split  in  two  subphases: 
determination  of  the  active  segments  (AS),  (i.e.,  those  signal 
portions  mainly  occupied  by  MUAPs),  and  fine  estimation  of 
the  span  of  these  segments. 

B.l.  Determination  of  AS.  Several  steps  are  carried  out: 

-  DWT  computation  of  the  EMG  signal  and  selection  of  the 
scale  containing  most  energy  (Fig.  2. a); 

-  First  splitting  of  the  signal  in  active  and  BL  segments:  as 
BL  samples  normally  dwell  in  a  narrower  amplitude  range 
than  MUAP  samples,  a  histogram  with  the  EMG  samples  (in 


the  wavelet  domain)  is  built  and  those  samples  whose 
amplitudes  appear  in  the  least  frequently  bin  are  assigned  to 
AS,  while  the  rest  are  assigned  to  BL  segments  (see 
thresholds  in  Fig.  2. a).  Histograms  of  20  to  100  bins  were 
tested,  finding  40-bin  histograms  as  the  most  convenient. 

Determination  of  maxima  and  minima:  maxima  and 
minima  in  AS  are  detected  and  considered  belonging  to 
different  MUAPs  or  superposition  of  MUAPs  when  a 
maximum  is  followed  by  a  minimum  or  vice  versa  and  when 
these  are  at  least  At  ms  apart.  We  experimentally  set  up  a 
suitable  relation  between  this  delay  and  the  EMG  activity: 
(low  activity:  At=9ms,  low-mid  activity:  At=llms,  mid-high 
activity:  At=13ms,  high  activity:  At=17ms). 

-  AS  extraction:  initial  and  final  AS  points  are  estimated  by: 

begining segment  =  firstt  -  f first  t  - end M)  (1) 

end segment  =  endi  + 1  {first  M  -  end , )  (2) 

where  firsts  and  endt  refer  respectively  to  the  first  maximum  or 
minimum  and  the  last  maximum  or  minimum  pertaining  to  a 
certain  MUAP  or  MUAP  superposition,  and  i  indexes  the 
different  MUAP  or  MUAP  superposition. 

B. 2.  Fine  estimation  of  the  initial  and  final  points  of  the 
identified  AS  (Fig.  2.b).  Several  steps  are  carried  out,  some  of 
them  equivalent  to  those  in  Subphase  B.l : 

-  DWT  computation  of  the  EMG  signal  and  selection  of  the 
scale  containing  most  energy. 

-  Second  splitting  of  the  signal  in  active  and  BL  segments:  a 
new  histogram-based  splitting  process,  similar  to  the  one 
previously  described  is  carried  out,  obtaining  more  restrictive 
AS.  Here  a  40-bin  histogram  was  also  chosen. 

-  Isolated  peak  elimination:  due  to  the  MUAP  morphology, 
the  maxima  and  minima  of  the  DWT  signal  should  appear 
alternately.  Peaks  altering  this  disposition  are  removed. 

-  Look  for  relevant  unconsidered  peaks:  some  low 
amplitude  MUAPs  may  have  skipped  from  the  previous 
processes.  They  are  search  for  at  either  side  of  the  AS,  using 
wavelet  domain  amplitude  thresholds,  peak-to-peak 
separation  time  and  sign  alternation  restrictions. 

Determination  of  the  potential  beginning:  the  onset  of  the 
first  maximum  or  minimum  in  the  wavelet  transform 
provides  the  beginning  of  the  MUAP  or  MUAP  superposition 
[7].  The  DWT  causes  an  artificial  delay  of  2j_1-2  samples 
which  has  to  be  discounted. 

-  Determination  of  the  potential  end:  the  end  point  of  the 
last  peak  is  taken  as  the  end  of  the  potential. 

C.  BL  analysis  and  spectral  characterization 

BL  spectral  characterizing  is  carried  out  in  this  phase,  by 
means  of  the  following  steps: 

-  Averaging  of  free  BL  segments:  The  detected  free-BL 
segments  may  still  contain  secondary  potentials  (low- 
amplitude,  smooth  potentials  coming  from  distant  MUs)  as 
well  as  high-frequency  noise  from  diverse  origins.  To  reduce 
the  influence  of  these  artifacts  in  BL  estimation,  consecutive 
samples  of  these  segments  are  averaged  (Fig.  2.c).  After 


testing  several  intervals  lengths  (between  3  and  20  ms), 
intervals  of  1 0  ms  were  chosen  as  the  most  convenient. 

-  Interpolation:  the  previously  averaged  points  are 
interpolated  by  means  of  cubic  splines,  resulting  a  curve  with 
the  appearance  of  a  true  BL  fluctuation  (Fig.  2.c). 

-  BL  spectral  characterization:  AR  spectral  estimation  [6] 
is  applied  to  the  interpolated  signal  giving  a  smooth  and  high 
resolution  power  spectral  density  (psd)  (Fig.  3).  As  expected, 
the  resulting  psd  corresponds  to  a  low  frequency  signal  (BL 
estimate).  Its  3-dB  bandwidth  is  then  obtained. 

D.  Filtering  (Fig.  2.d) 

For  the  final  BL  removal,  the  EMG  signal  passes  a  high- 
pass  filter  with  a  cut-off  frequency  equal  to  the  previous  3-dB 
bandwidth.  A  linear  phase  FIR  filter  is  used  so  to  preserve 
phase  relations  among  different  signal  components.  We  used 
Remez  interchange  algorithm  for  the  design  of  this  filter  [6]. 


IV.  Results 


A.  Merit  figures 

Quantitative  evaluation  is  needed  to  compare  BL  removal 
methods.  Two  merit  figures  (F  and  N)  are  devised  for 
measuring  the  degree  of  BL  fluctuation  in  repeated  discharges 
of  the  same  MU.  They  are  calculated  as  follows: 

1)  All  the  waveforms  in  the  EMG  record  corresponding 
to  non-overlapped  MUAPs  are  manually  selected. 

2)  These  discharges  are  given  time  origins  so  that 
correlation  among  them  is  maximum. 

3)  Let  Yk={yk(l ),  yk(2 ),  .  .,  yk(n) }  be  the  discharge  k  of 
the  set  of  n  discharges,  where  yk(t),  is  the  t  sample  of  Yk  . 
The  two  proposed  merit  figures  are  defined  as: 


F 


f 

=  sd 
k  V 


mean 

t 


A 


(3) 


(first  the  temporal  mean  of  every  discharge  is  calculated; 
then  the  standard  deviation  of  all  these  means  is  computed). 

N  =  mean  sd  (jq  (t\ ...,  ym(t)\  (4) 

t  \  k  ) 

(standard  deviation  across  different  discharges  is  calculated 
for  every  sample  time;  the  resulting  set  of  values  are  then 
averaged). 

BL  removal  methods  can  be  compared  using  these  figures: 
lower  F  and  N  values  are  attained  when  lower  BL  fluctuation 
remains,  indicating  thus  better  performance. 

To  first  test  the  validity  of  F  and  N ,  they  were  used  to 
measure  two  different  sets  of  simulated  signals.  The  first  set 
was  made  out  of  3  different  MUAP  trains  to  which  noise  (33 
dB  SNR  level)  and  an  artificial  BL  of  varying  amplitude 
(from  A  to  A/6)  and  low  frequency  range  (5Hz)  were  added. 
Table  I  shows  the  obtained  F  and  N.  Their  increasing  values 
with  BL  growing  amplitude,  credited  these  figures  as  good 
indicators  of  the  degree  of  BL  fluctuation.  In  the  second  set, 
the  same  MUAP  trains  were  contaminated  by  additive  noise 
(33  dB  SNR  level)  and  an  artificial  BL  of  different  frequency 
ranges  (5,  7,  10,  15  and  19  Hz).  Although  F  and  N  also 
increased  with  BL  frequency  content  (Table  I),  this  tendency 
was  not  so  notorious,  not  contradicting  the  previous  claim. 


Samples 

Figure  2. a.  DWT  and  amplitude  histogram-based  thresholds. 


Figure  2.b.  Active  segments  (separated  by  dashed  and  continuous  lines). 


Figure  2.c.  Free-BL  averages  (x)  and  interpolated  curve. 


Figure  3.  EMG  signal  spectrum  (a),  estimated  BL  (difference  between  initial 
signal  and  final  filtered  signal)  spectrum  (b),  AR  model  spectrum  (c),  and 
high-pass  filter  frequency  response  (d). 

B.  Results  from  simulation  and  real  signals 

Table  II  shows  the  mean  and  standard  deviation  of  F  and  N 
values  obtained  by  the  conventional  method  [1]  and  ours, 
when  applied  to  the  sets  of  real  and  simulated  signals.  The 
significant  lower  values  manifest  the  superior  performance  of 
our  method. 

C.  Visual  assessment. 

Visual  inspection  of  the  analyzed  signals  proved 
satisfactory  results  (compare  Fig.  1  and  Fig.  2.d),  except  in 
the  extreme  cases  of  high  activity  level  or  the  presence  of 
potentials  with  unusually  long  tails. 

IV.  Discussion 

It  can  be  noticed  that  the  differences  between  successive 
occurrences  of  the  same  bioelectric  phenomenon  (i.e.,  MUAP 
discharges)  are  lower  when  our  method  is  applied  than  when 
simply  removing  a  constant  BL  level.  From  this  we  infer  that 
BL  cause  an  artifactual  fluctuation  that  our  method  is  able  to 
counterbalance  in  some  extent. 

Table  i 

F  and  N  values  of  the  three  motor  units  (MU  1,  2  and  3)  composing  the 
simulated  EMG  records  with  varying  degrees  of  BL  amplitude  and  frequency 

content. 


MU  1 

MU  2 

MU  3 

F 

N 

F 

N 

F 

N 

BL 

amp. 

A/6 

40,5 

41,8 

11,9 

15,0 

26,9 

28,6 

A/5 

48,1 

49,2 

14,2 

17,1 

32,3 

33,7 

A/4 

59,4 

60,4 

17,8 

20,6 

40,4 

41,6 

A/3 

78,4 

79,12 

23,7 

26,6 

53,9 

55,0 

A/2 

116,3 

117,0 

35,5 

38,8 

80,9 

81,9 

A 

230,0 

230,7 

70,9 

76,3 

161,8 

163,2 

BL 

freq. 

5 

59,5 

60,4 

17,8 

20,6 

40,4 

41,6 

7 

51,0 

52,3 

36,1 

37,2 

76,5 

77,1 

10 

53,3 

55,7 

77,6 

79,8 

83,1 

84,3 

15 

74,9 

79,1 

54,9 

71,3 

99,2 

101,7 

TABLE  II 

F  and  N  mean  and  standard  deviation  values  obtained  by 
the  conventional  method  (Prig.)  and  ours  (New). 


F 

N 

Mean  (std) 

T-test 

Mean  (std) 

T-test 

Real 

Orig. 

143,7(260,1) 

P< 

387,5  (325,8) 

P< 

signals 

New 

129,9  (254,8) 

0,01 

381,9  (328,1) 

0,01 

Simul. 

Orig. 

59,8  (33,4) 

P< 

77,5  (54,8) 

P< 

Signals 

New 

37,0  (23,4) 

0,001 

53,76  (54,7) 

0,010 

When  the  EMG  signal  contains  waveforms  with  small-slope 
and  long  tails,  determination  of  their  ends  points  may  turn 
hazardous.  In  such  cases  averaging  cannot  fully  cope  with 
these  terminal  waveform  portions,  and  spectral 
characterization  becomes  inaccurate. 

In  the  case  of  high  EMG  activity  the  method  efficacy  is  also 
reduced:  as  this  activity  grows,  free-BL  segments  are  sparser 
and  shorter,  splines  curvature  is  abnormally  high  in  the 
interpolated  curve,  and  the  BL  course  is  not  followed  too 
precisely,  yielding  a  final  distorted  estimation  of  the  BL 
spectrum.  However,  this  limitation  is  not  too  problematic  as 
EMG  signals  with  such  activity  levels  are  unacceptable  for 
clinical  practice. 

The  proposed  BL  removal  method  can  be  used  for  obtaining 
more  precise  measurements  of  the  conventional  MUAP 
parameters  [1],  [2],  although  the  potential  benefit  of  this 
enhancement  has  still  to  be  explored.  Moreover,  being  a 
method  designed  to  be  applied  on  free-run  EMG  signals,  it 
may  eventually  be  implemented  as  part  of  EMG 
decompositions  systems  for  automatic  MUAP  extraction. 

V.  Conclusions 

The  proposed  method  for  BL  removal,  based  in  frequency 
characterization  of  MUAP-free  EMG  segments  has  proven 
superior  to  conventional  procedures,  that  consider  the  BL  to 
be  constant  along  the  MUAP  record  or  the  whole  EMG 
signal.  Its  main  limitation  is  the  presence  of  MUAPs  with  a 
long  tail  that  cannot  be  completely  put  aside  from  the  BL 
segments,  altering  thus  the  BL  frequency  estimation. 

Its  potential  use  for  enhancing  MUAP  parameter 
measurements  and  real-time  acquisitions  appears  promising. 
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